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On numerical study of the discrete spectrum of a 
two-dimensional Schrodinger operator with soliton 

potential 

A.N. Adilkhanov * I.A. Taimanov ^ 


Abstract 

The discrete spectra of certain two-dimensional Schrodinger oper¬ 
ators are numerically calculated. These operators are obtained by the 
Moutard transformation and have interesting spectral properties: their 
kernels are multi-dimensional and the deformations of potentials via 
the Novikov-Veselov equation (a two-dimensional generalization of the 
Korteweg-de Vries equation) lead to blowups. The calculations sup¬ 
ply the numerical evidence for some statements about the integrable 
systems related to a 2D Schrodinger operator. The numerical scheme 
is applicable to a general 2D Schrodinger operator with fast decaying 
potential. 


1 Introduction 

This paper concerns numerical calculation of the discrete spectra of some 
two-dimensional Schrodinger operators with soliton potential. These poten¬ 
tials were found in mm by using the Moutard transformation well-known 
in surface geometry and have interesting spectral properties. The numeri¬ 
cal scheme is based on the Galerkin method and can be applied to general 
operators with fast-decaying potentials. 

The Darboux transformation (the ladder method) allows us to construct 
the integrable one-dimensional Schrodinger operators 

(f , . 
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i.e. the operators whose spectra and eigenfunctions are explicitly described. 
To each integrable operator H and a solution ui of the equation Huj = 0 
there corresponds the new integrable operator H of the same form. Namely 
this correspondence sends each solution of the equation = Eip to a 
function (p = such that Hip = Ep, where A is some first order differential 
operator independent of ip. Note that this correspondence is invertible. 

For instance, the quantum harmonic oscillator is exactly the operator 
for which this transformation results in a shift of the potential u{x) by a 
constant. This observation was used by Dirac for finding the spectrum and 
eigenfunctions of the operator [3] . The successive application of the Darboux 
transformation, starting with the trivial potential u = 0 , yields the rational 
sohtons 

In fact, the Darboux transformation was originally derived as a one¬ 
dimensional reduction of the Moutard transformation which acts on two- 
dimensional Schrodinger operators as follows: 

(92 52 

H =-Add + u{x,y) =^+u{x,y). (2) 

Furthermore to each operator H and a solution uj of the equation Huj = 0 
there corresponds another operator H of this form, which gives an explicit 
procedure of constructing solutions to the equation Hp = ^ from solutions 
of the equation Hip = 0. The correspondence is one-to-one modulo In 

contrast to the Darboux transformation this method gives no information 
on solutions of the equation 

Hp = Ep 

with i? 7 ^ 0 . 

For the potentials derived by the Moutard transformation it is interesting 
to know something beyond the “zero energy level” E = Q. In this article we 
demonstrate how to achieve this aim numerically: 

• we calculate the discrete spectra of the Schrodinger operators with 
nontrivial kernel which were found in |T]; 

• we describe the dynamics of the discrete spectra of the operators whose 
potentials depend on t, the temporary variable, and give a blowing-up 
solution of the Novikov-Veselov equation [2]. 

We think that the numerical approach can be rather helpful for under¬ 
standing the spectral problems, will lead to formulating mathematical state¬ 
ments that are based on experiments and can be valuable for further proofs. 
Some of these conclusions are presented as final Remarks. 
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2 The Moutard transformation and some of its ap¬ 
plications 

2.1 The transformation 


Let H be of the form ([2]) and let w be a solution to the equation 

Hlo = (—A + u)ijJ = 0. 


Then the Moutard transformation of H is defined as 


H = — A + u — 2A logo; = —A — u + 2 


^x+^y 

■ 


Straightforward computations show that if ip satisfies the equation Hip = 
0 , then the function defined from the system 


, {(Jj(p)y = U? 

satisfies the equation 
If ip meets jS]), then 


(3) 


Hip = 0. 
C 


ip + 


UJ 


satisfies the same equation for every constant C. 

If the potential u = u{x) depends only on x and w = /(x)e'^^, then 

d? 




and the Moutard transformation reduces to the Darboux transformation of 
Hq defined by /: 


H = Ha - 


5y2 


H = Ha- 


dy^ 


Originally, the Moutard transformation was introduced for the hyperbolic 
operators 

+ u, 

and was used, for instance, for constructing new negatively curved surfaces 
in from those already available. However, since the procedure is formally 
analytical we may put P, = z and r] = z and derive its version for Schrodinger 
operators. 
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2.2 Two-dimensional Schrodinger operators with nontrivial 
kernel 

The Moutard transformation applied to the trivial potential u = 0 does not 
give a nonsingular fast decaying potential. However such a potential for 
which the scattering operator is well defined may be achieved by the double 
iteration of the Moutard transformation. That was done in |T] whose main 
results are as follows: 

• the potential (see Fig. [T]) 

5120(1 + 8x + 2?/ + 17x2 + 17y2) 

(160 + 4(x2 + y2)(i ^4x + y) + 17(x2 + y2)2)2 “ ^ 

5120|l + (4-i)z|2 
"(160 + |.2|2|2 + (4-z)z|2)2 

is nonsingular and decays as while the corresponding Schrodinger 



Figure 1: The potential u of the form(|T|) 
operator H of the form ([2]) has at least two-dimensional kernel (see Fig. 

ED: 

= H^2 = 0 

with 

^ X + 2x2 + xy — 2?/2 

^ 160 -b 4x2 -b 4y2 -b 16x3 -b 4x2?/ -b 16xy2 + 4y3 _|_ _)_ y2^2 ’ 

^ 2x + 2y + 3x2 _|_ iQxy — 2>y'^ 

^ 160 -b 4x2 -b 4y2 -b 16x3 -b 4x'^y + 16xy2 + 4y3 _|_ _)_ y2^2 ’ 

( 5 ) 
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(a) 


(b) ^-2 


Figure 2: The eigenfunctions Ti and T 2 of the form ([5]) 


• the potential 


where 


^ Fo{x,y) 
G{x,yf 


( 6 ) 


Fo{x,y) = -1280000(25+ 20x-287x2 + 60x^ + 1800x^-30?/-eOOxy- 

SOOx^y + 313?/2 + GOxy^ + SOOOx^y^ _ sOOy^ + I 8 OO/), 

G(x, y) = 40000 + lOOx^ + 40x^ - 387x'^ + 40x^ + 800x® - eOx^y- 

SOOx^y — 200x^?/ + lOOy^ + A0xy‘^ + 26x‘^y^ + SOx^y"^ + 2400 x ^?/2 _ 50 ?/^ — 

800x/ - 400x2/ + 413 / 4 Q^y 4 2400x2/ - 200y^ + 8 OO/, 

is nonsingular and decays as r“®, while the following functions ipi and 
?/’2 span a two-dimensional subspace in the kernel of the Schrodinger 
operator: 

H^i = H'^2 = 0 


for 


where 


Fiix,y) ^ F 2 {x,y) 

G{x,y)' ^ G{x,y) 


(7) 


Fi(x, y) = — lOx — 2x2 + 20x^ + 6xy + 60x^y + 2 / — 60x/ — 20/, 


F 2 {x-, y) = — lOx — 5x2 _ -Y{)y _|_ 2xy + 12Qx^y + 5 / — 40/. 
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Remark. Note that 

1) some good scattering theory is available for an A^-dimensional Schro- 
dinger operator if the potential of the operator meets the condition 


|m(x)| < 


const 


(1 + |x|)^+^’ 


X € 


t,N 


( 8 ) 


with e > 0 an arbitrary positive constant |1]; 

2) for the one-dimensional operators meeting the Faddeev condition 


/ + 00 
-OD 


|ri(x)|(l -|- \x\)dx < oo 


the discrete spectrum is finite and negative; 

3) for > 5 it is easy to construct a potential of the form u{x) = 

where G is a nonsingular positive function coinciding with the fundamen¬ 
tal solution of the Laplace equation G{x) = outside some compact 

set, such that u{x) meets ([8]) and the corresponding operator has nontrivial 
kernel. 

Potentials m and Q are the first examples of the two-dimensional 
Schrodinger operators meeting ([8]) and having a nontrivial kernel. Simi¬ 
lar examples in dimensions N = 3 and N = 4 have been unknown until 
recently. 

2.3 Blowing up solutions of the Novikov—Veselov equation 

The Novikov-Veselov (NV) equation [5j has the form 

Ut = d^U + B^U + 3d{VU) + 3d{VU) = 0, 

BV = BU, 

and can be derived as the compatibility condition of the system 

Hip = 0, 

dtp = {d^ + B^ + 3Vd + 3VB)p 
where BV = BU, dV = BU , and 

— 1 ?y 

H = dB + u = -A- - 

4 4 


(9) 


( 10 ) 


is a two-dimensional Schrodinger operator. 
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Its one-dimensional reduction corresponding to U = U(x) and U = V is 
the Korteweg-de Vries equation 

Ut = -^Uxxx + QUUx- 


The KdV equation admits the Lax representation 


dL 

dt 


[L,A], 


where L the one-dimensional Schrodinger operator with potential U. There¬ 
fore, the KdV evolution preserves the whole spectrum of L. In contrast to 
the KdV case, the NV equation has the Manakov triple representation 


dH 

dt 


[H,A] + BH 


and preserves only the zero-level spectrum. 

The Moutard transformation of H can be extended to some transfor¬ 
mation of solutions of the NV equation to another solutions of the same 
equation. By using the transformation, in [2] the first blowing up solution 
of the NV equation was constructed which is as follows: 


U = 2dd log $, V = 2d‘^ log 
<h(x, y, t) = 3(x^ -|- y^) -|- 4(x^ -|- y^) -|- 30 — 12t. 


( 11 ) 


For t < t* the potential U is nonsingular and decays as r~^. As t ap¬ 
proaches from below, the solution starts to oscillate in a bounded domain 
and at t = t* = yI the function <I> vanishes at a couple of points at which 
U has singularities: V is a rational function and for t = the denominator 
vanishes at these points, but the nominator vanishes at certain lines passing 
through the same points (see Fig. ED- 

We recall that I = ^(^2 U dx dy is the first integral of the NV equation, 
which restricts the possible types of singularities. 

Remarks. 1) For the KdV equation, i.e. the one-dimensional reduction 
of the NV equation, such a blowup is impossible. For instance, a solution 
of the KdV equation with smooth initial data Uo{x) = o(|x|) as |x| —oo is 
regular and unique for all times t |6]. 

2) In [7] it was noted that the more general function U = 2ddlog{a{x^ -|- 
y^)^ -|- b{x^ + y^) + c — 46t), with a, b and c real constants, meets the NV 
equation and gives a similar blowing up solution. 
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Figure 3: The singularities of U 


3) The modified Novikov-Veselov (mNV) equation, which is a simi¬ 
lar two-dimensional version of the modified Korteweg-de Vries equation, 
also admits blowups which are of a different nature. For this equation 
I = J U'^ dxdy is the first integral in contrast to the NV equation whose 
evolution preserves J U dx dy. Recently in there was constructed a 

solution to the mNV equation on the whole plane which decays rather fast 
at infinity and is real-analytic outside the only singular point x = y = t = 0, 
while the first integral I equals Svr for t ^ 0 and jumps to 27r at t = 0. 

3 A numerical method 

3.1 The basic scheme 

We will briefly explain the scheme in general. 

Take an orthonormal basis {ek},k = 1,2,..., for ^ 2 ( 1 ^™): 

(cj|cfc) — jjk = 1,..., 

where {u\v) is the standard inner product on L 2 (M™'). 

Denote by the orthogonal projection to the finite-dimensional sub¬ 
space 14 spanned by the first k vectors: 

Pn '■ L2 ^ Vn = Cei CSn- 

Given an operator H : L 2 (M™') —)• L 2 (M™’) we consider its approximation of 
the form 

Hn = PnH :Vn^Vn 
and find the spectrum of Pin- 

eW = {Fij’") < • • • < 
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The operator Hn is represented in the basis ei,..., by the matrix 


Ajk = {ej\Hnek) = {ej\Hek), 1 <j,k <n, 


which is the principal minor of the matrix A describing operator H in the 
basis {cfc}. 

Let H be the Schrodinger operator 


dy'^ ) 


+ u{x,y) 


with a fast decaying potential. The number of eigenvalues of Hn counted with 
multiplicities is equal to n and grows as n —>■ oo. Hence the whole spectrum 
cannot be considered as approximation of the discrete spectrum of H 
which can be even empty. 

However it is known that if the potential decays rather fast, i.e., 


u = o 



as |x| —>■ oo; 


then the Schrodinger operator has no positive eigenvalues m- Hence 

the sequence with the excluded positive elements not very close to 
zero, i.e., the truncated sequence, can be considered as approximation of the 
discrete spectrum of H as n ^ oo. 

We do not supply any strong mathematical evidence to this claim. But 
similar arguments are widely used in computational chemistry (quantum 
chemistry and molecular dynamics) and the results ht experiments. In this 
article we present the results of a similar numerical study of the discrete 
spectra of potentials from 12.21 and 12.31 


3.2 A numerical scheme based on the Hermite functions 

The Hermite functions 

(- 1 )" 


‘fk = 


, A: = 0,1, 2,. 


^J2^k\y/x dx^ 
form the orthonormal basis for L 2 {M.), i.e. 

(^ipj\(pk) = 6jk, 

and are the eigenfunctions of the quantum harmonic oscillator 

d^ 


dx^ 


+ X ] ipk — {2k + l)(pk- 
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Let us consider in L 2 (M^) the orthogonal basis formed by the products 

^j,k{x,y) = y^j{x)(Pk{y) 


of the Hermite functions of x and y. In this basis the Schrodinger operator 
is described by the matrix 

Since is the eigenfunction of the two-dimensional quantum harmonic 
oscillator 

/ ^2 ^2 \ 

+ x^ + y‘^)(fj^k = ( ~ + x'^ + y^j (fj^k = 2(j + k + l)(pj^k, 

we have 

{(fj^k\H(fl^rn) = (V^j,/c|(-A + 'u)(^/,r^) = ( 99 ^,/, | (u-+ 2 (/+ m+ 1 )= 

^ y “ 1 “ + I)i5j7(5/cm* 

In the sequel let LI(Ar) denote the operator 

^{N) = P(n)HP[m) 


where 

P(Ar) : L2(K^) span{vPj,fc,0 < j,k < N - 1] 

is the orthogonal projection of L 2 (]R^) to the A^^-dimensional subspace spanned 
by the two-dimensional Hermite functions (pj^^ with 0 < j, /c < iV — 1. 

3.3 The Gauss-Hermite quadrature for matrix elements 

For calculating the matrix elements we use the Gauss-Hermite 

quadrature (see m) : the integral over M of a function of the form 

f{x) = e~^^hs{x), 

with hs{x) a polynomial of degree s, is evaluated by the quadrature formula 

y °9 n n 

/ f{x)dx= / h{x)dx = '^Wihs{x) = '^LXif{xi), 
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where Xi are the roots of the Hermite polynomial of degree n while Wi, uji 
are the corresponding weights: 

2”-in!v^ ^2 _ 2”-in!0F 

This quadrature is exact for all s up to 2n — 1. If /(x) = g{x)e~^^^‘^ with 
g{x) G L^(M) for which the coefficients if the Hermite expansion 

of g satisfy \ck\ < ^(1 + then the quadrature error is bounded by 

0{n-^). 

In the two-dimensional case 


CXD CO 

If f{x,y)dxdy 


— CO —CO 


EE LOiUjf{xi,yj). 

i=l j=l 


For calculating the roots of the Hermite polynomials we use the algorithm 
that is proposed in |12) . 


4 Numerical results 

4.1 The Schrodinger operators with nontrivial kernel 

1. The potential written down in (jl]) is nonpositive and vanishes exactly at 
one point. Therefore, by the Rayleigh principle, the discrete spectrum of the 
corresponding Schrodinger operator H is nontrivial and the first eigenvalue 
equals to 

Ai = inf 
Ibl=i 

It appears that except the five eigenvalues Ai,A 2 ,A 3 ,Aoi and A02 the 
eigenvalues of ll{Ar) are greater than 0.1 for N = 16,20,25,32, and 50. The 
exceptional five eigenvalues split into the two groups: 

a) the negative eigenvalues Ai,A 2 ,A 3 that are separated from 0; 

b) the positive eigenvalues Aqi and A02 that converge to the zero from 
above as N grows. 

The results of computations are given in Table [1] 

Denote by W the span of the eigenfunctions Ti and 4'2 of the form ([5]) 
and by ipi, -02, V's, V'oi) '002 the eigenvectors of 11(7V) corresponding to Ai, 
-^2, A3, Aoi, Ao 2- Let ai stand for the angle between 0oi and W, i = 1,2. 
The calculations in Table [2] show that 0oi and 0o2 converge to functions in 
IT as grows. 

Therefore we arrive at the following: 
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A 

= 256 

iV2 = 400 

= 625 

iV2 = 1024 

iV2 = 2500 

Ai 

-1.80934 

-1.8093 

-1.80935 

-1.80936 

-1.80936 

A2 

-1.09104 

-1.09134 

-1.09158 

-1.09163 

-1.09163 

A3 

-1.01904 

-1.01929 

-1.01927 

-1.01927 

-1.01927 

Aoi 

0.0205608 

0.0129299 

0.00811659 

0.00487366 

0.00194779 

Ao2 

0.0320718 

0.0199062 

0.0125 

0.00745718 

0.0029613 


Table 1: The approximate eigenvalues of the Schrodinger operator with the 
soliton potential (jH) 


Numerical Conclusion 1. The discrete spectrum of the operator H with 
potential 0 consists of five eigenvalues (their approximations are given in 
Table 1) of which three are negative. The kernel of H is two-dimensional. 


cos a 

= 256 

= 400 

= 625 

= 1024 

= 2500 

cos ai 

cos Oi2 

0.991287 

0.991362 

0.992226 

0.992084 

0.993199 

0.992921 

0.994215 

0.993901 

0.995846 

0.995526 


Table 2: The cosines of the angles between the eigenvectors corresponding 
to Aoi and A02 and the linear span of Ti and T2 

The fast decay of the coefficients of expansions of if a in ihJ = 

0,N — 1) is demonstrated by the graphs of these coefficients parameter¬ 
ized by i,j and calculated for N = 32 (see Fig. H- The graphs of these 
eigenvectors are presented in Fig. [3 
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Figure 4: The expansions of eigenvectors 





Figure 5: The eigenvectors 


2. The analogous calculations for potential ([6]), which is also nonpositive, 
shows that the eigenvalues of that may correspond to the eigenvalues 
of H split into the two groups: 
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a) the negative eigenvalues Ai,..., A5 separated from zero; 

b) the eigenvalues Aoi,Ao 2 converging to zero as N grows. 
The calculations of eigenvalues are presented in Table |3j 


A 

= 256 

N‘^ = 400 

= 625 

= 1024 

N‘^ = 2500 

Ai 

-2.83707 

-2.8501 

-2.84788 

-2.84807 

-2.84812 

A2 

-2.49391 

-2.51335 

-2.50945 

-2.50925 

-2.50919 

A3 

-2.23477 

-2.2435 

-2.24252 

-2.24291 

-2.24302 

A4 

-1.33696 

-1.34495 

-1.33922 

-1.33888 

-1.33882 

A5 

-1.31464 

-1.33525 

-1.33505 

-1.33493 

-1.33487 

Aoi 

0.008239 

0.000002 

0.000656 

0.000626 

0.000237072 

Ao2 

0.019057 

0.001244 

0.001943 

0.000966 

0.000248573 


Table 3: The approximate eigenvalues of the Schrodinger operator with the 
soliton potential Q 

The eigenvectors ipQi and ■002 corresponding to Aqi and A02 have to ap¬ 
proximate functions from the kernel. The cosines of the angles ai,i = 1,2, 
between the eigenvectors corresponding to Aqi and A02 and the linear span 
of Ti and '£'2 are presented in Tabled 


cos a 

= 256 

= 400 

iV2 = 625 

iV2 = 1024 

= 2500 

cos ai 

cos Oi2 

0.998691 

0.998792 

0.998852 

0.998910 

0.999245 

0.999273 

0.999499 

0.999505 

0.999764 

0.999765 


Table 4: The cosines of the angles between the eigenvectors corresponding 
to Aqi and A02 and the linear span of Ti and 'I'2 

Therefore, we arrived at 

Numerical Conclusion 2. The discrete spectrum of the operator H with 
potential m consists of seven eigenvalues (their approximations are given in 
Tahle\^ of which five are negative. The kernel of H is two-dimensional. 

4.2 The blowing up solution of the Novikov—Veselov equa¬ 
tion 

For the operator 

H = -A-AU 
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where U is of the form (fTT]) . i.e., a blowing up solution of the NV equation, 
we numerically find the five eigenvalues: 


Ai, A2, A3, Aoi, Ao2- 

At t = 0 the three eigenvalues Ai,A 2 ,A 3 are negative and lie rather far 
from zero, whereas Aqi and A02 are close to zero and approximate the zero 
eigenvalue. Indeed, the kernel is at least two-dimensional and includes the 
subspace that is spanned by 

2x‘^ — 2y‘^ + 2x + 2y 

“ 3x4 + 43.3 + e^2y2 3y4 + 4y3 + SQ - I2t’ 

“4xy 

^ ^ 3x4 _j_ 42-3 _|_ _j_ 3^4 _|_ 4^3 _)_ 30 — l2t 

Since we have an explicit description of the evolution of the potential we 
can calculate the dynamics of the eigenvalues, which is presented in Fig. [6j 



Figure 6: The approximate dynamics of the eigenvalues for N = 50 

The vertical line denotes the critical time T* = ^ at which the potential 
becomes singular. 

We see the following: 

• Up to t ~ 2.2 the two maximal eigenvalues stay closed to the zero 
as it has to be because the zero-level spectrum is preserved by this 
evolution. Therewith the negative spectrum evolves. 
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• After t ~ 2.2 the approximations to the zero eigenvalue substantially 
decrease going away from zero and therefore the calculations become 
unreliable. The correct description of the evolution of the discrete 
spectrum for t > 2.2 needs the Galerkin approximations with respect 
to higher-dimensional subspaces or a more precise numerical scheme. 


5 Final remarks 

The calculations were done by using the LAPACK package. 

The main points for us consist in 

1 ) demonstrating the possibilities of the method and we see, for instance, 
that the kernel is detected rather precisely; 

2 ) studying by numerical methods the spectral properties of some two- 
dimensional differential operators. 

Usually the results of calculations are compared with an integrable case. 
Since the lack of examples of fast-decaying two-dimensional potentials with 
explicitly described discrete spectrum we use for such a comparison the zero- 
level spectrum explicitly established for the potentials that are considered in 

IB. 

These calculations provide the numerical evidence for the following state¬ 
ments: 

1. The negative discrete spectrum is not preserved by the Novikov-Veselov 
Bow and it seems that all conservation laws are given by the spectral 
data of H which correspond to the zero energy level (see m- 

2. The Darboux transformation can be interpreted as a transformation 
that adds or removes a point from the discrete spectrum of a one¬ 
dimensional operator ©• The operators from 2.2 and 2.3 are derived 
by the double iteration of the Moutard transformation starting at the 
zero potential u = 0. The analogy between these transformations 
and the numerical calculations leads us to the conjecture that if the 
Moutard transformation respects the function classes of potentials and 
eigenfunctions then 

the dimension of the kernel of H is changed by one by the Moutard 
transformation and there are no such estimates related to the negative 
part of the discrete spectrum. 
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in particular, for the suggestion to apply the Galerkin method. 
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